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Abstract 

We investigate constrained optimal control problems for linear stochastic dynamical systems evolving in discrete time. We 
consider minimization of an expected value cost over a finite horizon. Hard constraints are introduced first, and then reformu- 
lated in terms of probabilistic constraints. It is shown that, for a suitable parametrization of the control policy, a wide class 
of the resulting optimization problems are convex, or admit reasonable convex approximations. 
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1 Introduction 

This work stems from the attempt to address the optimal infinite-horizon constrained control of discrete-time stochas- 
tic processes by a model predictive control strategy [1,14,12,13,15,29,2,4,7]. We focus on linear dynamical systems 
driven by stochastic noise and a control input, and consider the problem of finding a control policy that minimizes 
an expected cost function while simultaneously fulfilling constraints on the control input and on the state evolu- 
tion. In general, no control policy exists that guarantees satisfaction of deterministic (hard) constraints over the 
whole infinite horizon. One way to cope with this issue is to relax the constraints in terms of probabilistic (soft) 
constraints [25,26]. This amounts to requiring that constraints will not be violated with sufficiently large probability 
or, alternatively, that an expected reward for the fulfillment of the constraints is kept sufficiently large. 

Two considerations lead to the reformulation of an infinite horizon problem in terms of subproblcms of finite horizon 
length. First, given any bounded set (e.g. a safe set), the state of a linear stochastic dynamical system is guaranteed 
to exit the set at some time in the future with probability one whatever the control policy. Therefore, soft con- 
straints may turn the original (infeasible) hard-constrained optimization problem into a feasible problem only if the 
horizon length is finite. Second, even if the constraints are reformulated so that an admissible infinite-horizon policy 
exists, the computation of such a policy is generally intractable. The aim of this note is to show that, for certain 
paramctcrizations of the policy space [21,6,17] and the constraints, the resulting finite horizon optimization problem 
is tractable. 



An approach to infinite horizon constrained control problems that has proved successful in many applications is 
model predictive control [22]. In model predictive control, at every time t, a finite- horizon approximation of the 
infinite-horizon problem is solved but only the first control of the resulting policy is implemented. At the next time 
t + 1, a measurement of the state is taken, a new finite-horizon problem is formulated, the control policy is updated, 
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and the process is repeated in a receding horizon fashion. Under time-invariance assumptions, the finite-horizon 
optimal control problem is the same at all times, giving rise to a stationary optimal control policy that can be 
computed offline. 

Motivated by the previous considerations, here we study the convexity of certain stochastic finite-horizon control 
problems with soft constraints. Convexity is central for the fast computation of the solution by way of numerical 
procedures, hence convex formulations [8] or convex approximations [23,9] of the stochastic control problems are com- 
monly sought. However, for many of the classes of problems considered here, tight convex approximations are usually 
difficult to derive. One may argue that non-convex problems can be tackled by randomized algorithms [28,27,30]. 
However, randomized solutions are typically time-consuming and can only provide probabilistic guarantees. In par- 
ticular, this is critical in the case where the system dynamics or the problem constraints are time-varying, since in 
that case optimization must be performed in real-time. 

Here we provide conditions for the convexity of chance constrained stochastic optimal control problems. We derive and 
compare several explicit convex approximations of chance constraints for Gaussian noise processes and for polytopic 
and ellipsoidal constraint functions. Finally, we establish conditions for the convexity of a class of expectation-type 
constrains that includes standard integrated chance constraints [19,20] as a special case. For integrated chance 
constrains on Gaussian processes with polytopic constraint functions, an explicit formulation of the optimization 
problem is also derived. 

The optimal constrained control problem we concentrate on is formulated in Section 2. A convenient parametrization 
of the control policies and the convexity of the objective function are discussed at this stage. Next, two probabilistic 
formulations of the constraints and conditions for the convexity of the space of admissible control policies are 
discussed: Section 3 is dedicated to chance constraints, while Section 4 is dedicated to integrated chance constraints. 
In Section 5, numerical simulations are reported to illustrate and discuss the results of the paper. 



2 Problem statement 



Let N = {1, 2, . . .} and N = N U {0}. Consider the following dynamical model: for t £ N , 

x(t + 1) = Ax(t) + Bu(t) + w(t) 



(1) 



where x(t) G W 1 is the state, u(t) € K m is the control input, A e K nX ™, B e W ixm , and w(t) is a stochastic noise 
input defined on an underlying probability space (CI, P). No assumption on the probability distribution of the 
process w is made at this stage. We assume that at any time t € No, x(t) is observed exactly and that, for given 
x £ K", x(0) = x Q . 

Fix a horizon length JVeN. The evolution of the system from t = through t = N can be described in compact 
form as follows: 

x = Ax + Bu + Dw, (2) 

where 
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Let V : R( Ar + 1 )" xWm ^ K and 77 : R(JV+i)nxJVm _> with r e ^ be measurable functions. We are interested in 
constrained optimization problems of the following kind: 



inf E[V(x,u)\ 
subject to (2) and r](x,u) < 



(3) 



where the expectation E[-] is defined in terms of the underlying probability space P), % is a class of causal 

deterministic state-feedback control policies and the inequality in (3) is interpreted componentwise. 

Example 1 In the (unconstrained) linear stochastic control problem [3] , w is Gaussian white noise and the aim is 
to minimize 

"JV-l 

E (x T {t)Q(t)x{t)+u T {t)R(t)u{t)) +x T (N)Q(N)x{N) , 
. t=o 

where the matrices Q(t) G M™ x " and R(t) G flj mx ™ arc positive definite for all t, with respect to causal feedback 
policies subject to the system dynamics (1). This problem fits easily in our framework; it suffices to define 



V(x,u) 



-T -T 
X U 



M 



(4) 



with M = diag(Q(0), Q(l), . . . , Q(N),R(0), i?(l), . . . , R(N - 1)) > (the notation M > indicates that M is a 
positive definite matrix). In our framework, though, the input noise sequence may have an arbitrary correlation 
structure, the cost function may be non-quadratic and, most importantly, constraints may be present. 

Standard constraints on the state and the input are also formulated easily. For instance, sequential ellipsoidal 
constraints of the type 



x T (t) u T (t) 



S(t) 



x(t) 
«(*) 



< 1, 



t = 0,l,...,N-l, 



x T (N)S(N)x(N) < 1, 



with < S(t) G R(™+m)x(r l+ m) for f = 0, 1, . . . , JV - 1 and < S(N) G R nx ", are captured by the definition 

1 T 



rj (x,u) r)i(x,u) ■ ■ ■ rj N (x) 



where, for t = 0, 1, 
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and each matrix E! t is immediately constructed in terms of the S(t). Our framework additionally allows for cross- 
constraints between states and inputs at different times. 

2.1 Feedback from the noise input 



By the hypothesis that the state is observed without error, one may reconstruct the noise sequence from the sequence 
of observed states and inputs by the formula 



w(t) = x(t + 1) - Ax(t) - Bu(t), teN . 
In light of this, and following [21,17,6], we shall consider policies of the form: 

t-i 



(5) 



(6) 



i=0 
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where the feedback gains Gt,% G R mx ™ and the affine terms dt G K m must be chosen based on the control objective. 
With this definition, the value of u at time t depends on the values of w up to time t — 1. Using (5) we see that u(t) is 
a function of the observed states up to time t. It was shown in [17] that there exists a (nonlinear) bijection between 
control policies in the form (6) and the class of affine state feedback policies. That is, provided one is interested in 
affine state feedback policies, parametrization (5) constitutes no loss of generality. Of course, this choice is generally 
suboptimal, since there is no reason to expect that the optimal policy is affine, but it will ensure the tractability of 
a large class of optimal control problems. In compact notation, the control sequence up to time N — 1 is given by 

u = Gw + d, (7) 

where 
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(8) 



note the lower triangular structure of G that enforces causality. The resulting closed-loop system dynamics can be 
written compactly as the equality constraint 

x = Ax Q + B(Gw + d) + Dw. (9) 



Let us denote the parameters of the control policy by = (G, d) and write (xg, ug) to emphasize the dependence of 
x and u on 0. From now on we will consider the optimization problem 

inf nV(xg,ug)} (10) 
gee 

subject to (7), (9) and (11) 

ri(xe,u e )<0, (12) 

where is the linear space of optimization parameters in the form (8). 

Remark 2 With the above parametrization of the control policy, both ug and xg are affine functions of the param- 
eters (for fixed w) and of the process noise w (for fixed 0). Most of the results developed below rely essentially on 
this property. It was noticed in [5] that a parametric causal feedback control policy with the same property can be 
easily defined based on indirect observations of the state, provided the measurement model is linear. The method 
enables one to extend the results of this paper to the case of linear output feedback. For the sake of conciseness, this 
extension will not be pursued here. 

2.2 Optimal control problem with relaxed constraints 

In general, no control policy can ensure that the constraint (12) is satisfied for all outcomes of the stochastic input w. 
In the standard LQG setting, for instance, any nontrivial constraint on the system state would be violated with 
nonzero probability. We therefore consider relaxed formulations of the constrained optimization problem (10)— (12) 
of the form 

inf E[V(x e ,ug)] (13) 
gee 

subject to (7), (9) and (14) 
E[4>on(x0,He)]<O, (15) 

where <p : M. r -> R R , with ReN, is a convenient measurable function and the inequality is again interpreted compo- 
nentwise. For appropriate choices of </>, this formulation embraces most common probabilistic constraint relaxations, 
including chance constraints (see e.g. [23]), integrated chance constraints [19,20], and expectation constraints (see 
e.g. [26]). 

We are interested in the convexity of the optimization problem (13)— (15). First we establish a general convexity 
result. 
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Proposition 3 Let (Q,^, P) be a probability space, be a convex subset of a vector space and 2> C R be convex. 
Let 7 : 17 x — > £F and </? : 21 — > M 6e measurable functions and define 

J(e)=E[(po-f(u},e)). 

Assume that: 

(i) the mapping "/(uj, •) : — > R is convex for almost all ui G 17; 

(ii) </5 is monotone nondecreasing and convex; 
(hi) J(6>) is /mi£e /or all eQ. 

Then the mapping J : — + R is convex. 

Proof: Fix a generic w £ 17. Since 7(0;, 0) is convex in 9 and 93 is monotone nondecreasing, for any 9,9' E and any 
aG [0,1], 

p(7(w, a0 + (1 - a)6')) < tp(cry(w, 9) + (1 - a) 7 (w, 9')) . 

Moreover, since <p is convex, 

<p( ai {u, 9) + (1 - a) 7 (w, 0')) < a*>(7(w, »)) + (1 - "Mt^, 0')) • 
Since these inequalities hold for almost all cj G 17, it follows that 

E[p(7(w, a0 + (1 - a)0'))] < E[a^( 7 (a;, 0)) + (1 - a)ip(y(u, 9'))} 

= aE[ V ( 7 ( W> 0))] + (1 - a)E[p( 7 (« > , ))] J 

which proves the assertion. □ 

Assumption (iii) can be replaced by cither of the following: 

(hi') J(9) G R U {+00}, WO G 0. 
(iii") J (9) G R U {—00}, V0 G 0. 

Let us now make the following standing assumption. 

Assumption 1 V(x,u) is a convex function of (x,u) and E[V(xe, ug)] is finite for all 9 G 0. 
Proposition 4 Under Assumption 1, E[V(xg, ug)] is a convex function of 9. 

Proof: First, note that the set of admissible parameters 9 is a linear space. Let us write w(tu) xg{u>) and Ug(u>) to 
express the dependence of w, xg and ug on the random event w G 17. Fix w arbitrarily. Since the mapping 



Xg{uj) 




Ax + B(Gw(to) ^ 


- d) + Dw(lu) 


Ug(ui)_ 




Gw{u>) - 


hd 



is affine and the mapping (x, u) 1— > V{x, u) is assumed convex, their combination 9 1— ► V{xg, ug) is a convex function 
of 9. Then, the result follows from Proposition 3 with 7(0;, 9) = V(x$(ui), ug(u>)) and tp equal to the identity map. □ 

By virtue of the alternative assumptions (iii') and (hi") of Proposition 3, the requirement that E[V(xe,«e)] be hnite 
for all 9 may be relaxed. A sufficient requirement is that there exist no two values 9 and 9' such that J Xo (9) = +00 
and J Xa {9') = —oo. In particular, the result applies to quadratic cost functions of the type (4), with M > 0. 

In general, the relaxed constraint (15) is nonconvex even if the components rji : '^( N + 1 ) nxNm _> R ; with i = 
1, . . . ,r, of the vector function 77 are convex. In the next sections we will study the convexity and provide convex 
approximations of (15) for different approaches to probabilistic relaxation of hard constraints, i.e. for different choices 
of the function </>. 
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3 Chance Constraints 



For a given a £ ]0, 1[, we relax the hard constraint rj(xe,ug) < by requiring that it be satisfied with probability 
1 — a. Hence we address the optimization problem 

inf E[V(x e ,u e )} (16) 
eee 

subject to (7), (9) and (17) 

F(rj(xe,u e ) < 0) > 1 -a. (18) 

The smaller a, the better the approximation of the hard constraint (12) at the expense of a more constrained 
optimization problem. This problem is obtained as a special case of Problem (13)-(15) by setting R = 1 and defining 

4> as 

r 

(t>cc(ri) = 1 - JJlj-oo.o]^) - a, 

i=l 

where li_oo,o](") is the standard indicator function. We now study the convexity of (18) with respect to 9. 
3.1 The fixed feedback case 

First assume that the feedback term G in (7) is fixed and consider the convexity of the optimization problem (16)- 
(18) with respect to the open loop control action d. That is. for a given G in the form (8), the parameter space 
becomes the set {(G, d), Vd € R Nrn }. For the given G and i = 1 , . . . , r define 

gi(d, w) = r)i(Ax + B(Gw + d) + Dw, Gw + d), 

where ?/; : ]g( Ar + 1 )™+ Ar " 1 _ > ]£ i s the i-th element of the constraint function r\. Define 

p(d) = F[ gi (d, w)<0,... g r (d 7 w) < 0] (19) 

and &cc — {d ■ p(d) > 1 — a}. Observe that &cc corresponds to the constraint set dictated by (17)-(18) when G 
is fixed. 

Proposition 5 Assume that w has a continuous distribution with log-concave probability density and that, for i = 
1, ...,r, gi : ^ Nm + Nn _ > j s quasi-convex. Then, for any value of a €] 0, 1[, &cc * s convex. As a consequence, 
under Assumption 1 and for any a £] 0, 1[, the optimization problem 

inf E[V(xe,ue)] 
eee 

subject to (7), (9) and p(d) > 1 — a 

is convex. 

Proof: It follows from [24, Theorem 10.2.1] that (19) is a log-concave function of d, i.e. the mapping d i— » logp(d) is 
concave. Since log is a monotone increasing function, we may write that ^cc = {d '■ log p(d) > log(l — a)}. Hence, 
is a convex set. The convexity of the optimization problem follows readily from Assumption 1. □ 



Among others, Gaussian, exponential and uniform distributions are continuous with log-concave probability density. 
As for the functions gi, one case of interest where the assumptions of Proposition 5 are fulfilled is when n(x,u) is 
affiiie in x and u. This is the case of polytopic constraints, which will be treated extensively in the next section. 
Apparently, this convexity result cannot be applied to the ellipsoidal constraints treated subsequently in Section 3.3, 
nor can it be extended to the general constraint (18) with both d and G varying. Loosely speaking, the latter is 
because the functions gi are not simultaneously quasi-concave in w and G. In the next sections we will develop 
convex conservative approximations of (18) for various definitions of n. 
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3.2 Polytopic Constraint Functions 



Throughout the rest of Section 3 we shall rely on the following assumption. 

Assumption 2 w is a Gaussian random vector with mean zero and covariance matrix S > 0, denoted by w ~ 
Af(0,t). 

Polytopic constraint functions 

rj(x e , us) = T x x e + T u u e - y, (20) 

where T x £ K rx ( A '+ 1 )™ ; T u e R rxArm i and y e R r , describe one of the most common types of constraints. In light 
of (7) and (9), 

T]{xg,Ug) = hg + PgW, (21) 

where he = {T x Ax -y) + (T X B + T u )d, and Pg = (T X D + (T X B + T U )G). It is thus apparent that r](x e , u e ) is affine 
in the parameters 9. Yet, in general, constraint (18) is nonconvex. We now describe three approaches to approximate 
constraints (18) that lead to convex conservative constraints. 

3.2.1 Approximation via constraint separation 

Constraint (18) requires us to satisfy rj(xg,ug) < with probability of at least 1 — a. Here rj £ W and the inequality 
77 < is interpreted componentwise. One idea is to select coefficients on £ (0, 1) such that X3»=i a i = a an( i ^° satisfy 
the inequalities ¥(r/i < 0) > 1 — on, with i = 1, . . . , r. Note that this choice is obtained in (15) by setting R — r and 
4>i(yi) — 1 — 1]— 00,0] [Wi) ~~ a i where, for i = 1, . . . ,r, <f>i : W — > K denotes the i-th component of function </>. 

Let /ij^ and PJq be the i-th entry of hg and the i-th row of Pg, respectively, and let be a symmetric real matrix 
square root of £. 

Proposition 6 Let G]0, 1[. Under Assumption 2, the constraint 

F(r]i(xe,ue) < 0) > 1 - a h 
with rj defined as in (20), is equivalent to the second-order cone constraint in the parameters 9 £ 



fa 



< 



where /3j = a/2 erf 1 (1 — 2aj) and erf 1 (-) is the inverse of the standard error function erf(a;) = e u du. As 

a consequence, under Assumption 1 and if a\ + . . . + a r = a, the problem 



inf E[V(xg,ug)} 
see 

subject to (7), (9) and h^g + 
is a convex conservative approximation of Problem (16) -(18) 



£3R 



<0, i = l,. 



Proof: From Eq. (21) we can write r\i = hi t g + PjgW. Since w ~ 7V(0, X), 77.; is also Gaussian with distribution 
N{hi t e, P/gSPj.e). It is easily seen that, for any scalar Gaussian random variable X with distribution Af(n, <r 2 ), 



V{X < 0) > 1 - a <^=> /i + (3a < 0, 



< 0, 



where (3 = y/2 erf 1 (1 — 2a). Hence the constraint ¥(r]i(xg, ug) < 0) > 1 — ai is equivalent to hi.g + j3i 

where (3i = v / 2erf 1 (l — 2a,). Since hi_g and P^g are both affine in the parameters 9 = (G, J), the above constraint 
is a second-order cone constraint in 9. It is easy to see that this choice guarantees that P(n(xg, ug) < 0) > 1 — a. □ 



7 



3.2.2 Approximation via confidence ellipsoids 



The approach of Section 3.2.1 may be too conservative since the probability of a union of events is approximated by 
the sum of the probabilities of the individual events. One can also calculate a conservative approximation of the union 
at once. Constraint (18) with n(xg,ug) as in (21) restricts the choice of Pg and hg to be such that, with a probability 
of 1 — a or more, the realization of random vector rj(xg, ug) lies in the negative orthant rj(xg, ug) < 0. In general, it 
is difficult to describe this constraint explicitly since it involves the integration of the probability density of r\ over 
the negative orthant. However, an explicit approximation of the constraint can be computed by ensuring that the 
100(1 — a)% confidence ellipsoid of i] is contained in the negative orthant. Fulfilling this requirement automatically 
implies that the probability of rj(xg,ug) < is strictly greater than 1 — a. 

Since w ~ Af(0, E), it follows that rj(xg,ug) = hg + Pgw ~ Af(hg, Eg), with Eg := PgY.P e T . Consider the case where 
Ug is invertible. Define the r-dimcnsional ellipsoid 

£(hg,£g,(3) = { V eR r : | fa - hgftg 1 (n - hg) < (3 2 ) , (22) 

where (3 > is a parameter specifying the size of the ellipsoid. Notice that, in general, Eg is invertible when r < Nn 
(i.e. the number of constraints is less than the total dimension of the process noise). If r > Nn, as there arc Nn 
independent random variables in the optimization problem, the following result still holds with Nn in place of r. 

Proposition 7 Let a € ]0, 1[. Under Assumption 2, the constraint F[n(xg, ug) < 0] > 1 — a with r/ defined as in (20) 
is conservatively approximated by the constraint 

£(he,te,f3(a))(l(-^,0] r , (23) 

where (3(a) = y/F- 1 ^ — a) and F(-) is the probability distribution function of a x 2 random variable with r degrees 
of freedom. Moreover, (23) can be reformulated as the set of second- order cone constraints 

h il o+/9(a)||E'Pi, fl | <0, i = l,...,r. (24) 

As a consequence, under Assumption 1, the problem 

inf E[V(xg,ug)} 
subject to (7), (9) and (24) 

is a convex conservative approximation of Problem (16)— (18). 

Proof: Since r](xg,ug) ~ M(hg, £#), the random variable 

(n(xg,Ug) - hg) T '£,~ 1 (ri(xg,Ug) - hg) 

is x 2 with r degrees of freedom. Then, choosing j3 such that F((3 2 ) = 1 — a guarantees that £{hg, Eg, /3(a)) is the 
100(1 — a)% confidence ellipsoid for n(xg, ug). Finally, under (23), P[r](xg, ug) < 0] > P[n(xg, ug) <G £(hg,Y,g, (3(a))] = 
1 — a, which proves the first claim. To prove the second claim, note that (22) can alternatively be represented as 

£(hg,Eg,(3) ={r i £R r \l ] = hg + AIgU, \\u\\ <l}, (25) 

_ 1_ 

where M$ = /3E| . Since r\ < if and only if ejr) < 0, i = 1, ■ • ■ , r, where the denote the standard basis vectors 
in R r , we may rewrite (23) as 

ef(hg + M e u) < V||u||<l, 

or equivalently 

sup e[(hg + M 9 u) < 
IMKi 
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Fig. 1. Values of j3 as a function of the horizon length N for the ellipsoidal method (dashed line) and the constraint separation 
method (solid lines). It is assumed that r = N -r- (n + m), where n = 5 is the dimension of the system state and of the process 
noise, m — 2 is the dimension of the input u. r reflects the number of constraints per stage; we take f = 2, 10, 100, 1000. The 
value of (3 is invariant with respect to f (single dashed line) for the ellipsoidal method, while it increases with f (multiple solid 
lines) for the constraint separation method. 

for i =!,■■■ ,r. For each i, the supremum is attained with u = Mje^/ M^ei ; therefore, the above is equivalent 

_ i 

to efhg + Mg'ej < 0. Clearly, efhg = hi t $. Moreover, since Mg = , we have 



Therefore, constraint (23) reduces to (24). Since the variables h^g and P^g are affine in the original parameters 
8 = (G, d), this is an intersection of second order cone constraints. As a result, under the additional Assumption 1, 
the optimization of E[V(xg,ug)} with (24) in place of (18) is a convex conservative approximation of (16)-(18). □ 



Mj ei \\ = pJefZge, = pJeJPg^Pje, = [3 



3.2.3 Comparison of constraint separation and confidence ellipsoid methods 
In light of Propositions 6 and 7, both approaches lead to formulating constraints of the form 

Hi + (3cJi < 0, 

with i = 1, . . . , r, where [ii an d &i are the mean and the standard deviation of the scalar random variables r)i(xg, ug), 
and smaller values of the constant (3 correspond to less conservative approximations of the original chance constraint. 
For a given value of a, the value of (3 depends on the number of constraints r in a way that differs in the two cases. In 
the confidence ellipsoid method, in particular, the value of (3 is determined by Nn (the total dimension of the process 
noise) when r > Nn. In Figure 1, we compare the growth of (3 in the two approaches under the assumption that r 
grows linearly with the horizon length N. (For the constraint separation method we choose a\ = . . . = a r = a/r, 
so that the value of (3 is the same for all constraints.) The increase of (3 is quite rapid in the confidence ellipsoid 
method, which is only effective for a small number of constraints. An explanation of this phenomenon is provided 
by the following fact, that is better known by the name of (classical) "concentration of measure" inequalities; proofs 
may be found in, e.g., [10]. 
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Proposition 8 Let r^iy be the r -dimensional Gaussian measure with mean h and (nonsingular) covariance 
i.e., 

Then for ee ]0,1[, 



w r h , s ,({e| ne -ftii E /-i < v^(i-e)}) <c-^. 



The above proposition states that as the dimension r of the Gaussian measure increases, its mass concentrates in an 
ellipsoidal shell of 'mean-size' yfr. It readily follows that since r)(xg, ug) is a r-dimensional Gaussian random vector, 
its mass concentrates around a shell of size \fr. Note that the bounds corresponding to (i) and (ii) of Proposition 8 in 
the case of r](xg,ug) are independent of the optimization parameters 9; of course the relative sizes of the confidence 
ellipsoids change with 9 (because the mean and the covariance of rj(xg,ug) depend on 6), but Proposition 8 shows 
that the size of the confidence ellipsoids grow quite rapidly with the dimension of the noise and the length of the 
optimization horizon. Intuitively one would expect the ellipsoidal constraint approximation method to be more 
effective than the cruder approximation by constraint separation. Figure 1 and Proposition 8 however suggest that 
this is not the case in general; for large numbers of constraints (e.g. longer MPC prediction horizon) the constraint 
separation method is the less conservative. 

3.2-4 Approximation via expectations 

For any r-dimensional random vector 77, we have 1 — F(r] < 0) = E [l — Y\i=i l]-oo,o]( 7 7i)] ■ Using this fact one 
can arrive at conservative convex approximations of the chance-constraint (18) by replacing the function in the 
expectation with appropriate approximating functions. For ti > 0, i = 1, . . . , r, consider 

r 

<P(v) =^2exp(Urn). 



Lemma 9 For any r-dimensional random vector rj, £[93(77)] > 1 — P[t7 < 0]. 

Proof: For every fixed value of 77, it holds that 93(77) > l — J^^i l]-oo,o]( 7 7i)- Hence £[93(77)] > E[l — 11*= 1 l]-oo,o](*7t)] = 
1-P[77<0]. □ 



Proposition 10 Under Assumption 2, for 77 defined as in (20), it holds that 

r / t 2 
W.[ip(r](xg,Ug))] = ^]cxp (tih it g + ^-||S'Pi,« 



As a consequence, under Assumption 1 and for any choice of ti > 0, i = 1, . . . ,r, the problem 



inf E[V(xg,u e )] 
gee 



subject to (7), (9) and ^exp (t^.g + -^-||S^P 4! e|| 2 ^ < a 

i=l 

is a convex conservative approximation of Problem (16)— (18). 
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Proof: It is easily seen that, for any r-dimensional Gaussian random vector 77 with mean /i and covariance matrix E, 
and any vector c G M. r , E[exp(c T 77)] = cxp(c T /t + ic T £'c). Let us now write 7/0 in place of r)(xg,ug) for shortness. 
By the hypotheses on w, in the light of (21), 770 is Gaussian with mean hg and covariance Pg'SPj '. Then, for a vector 
Q G M r with zero entries except for a coefficient ij in the z-th position, 



E[cxp(cf 7/e )] = exp (cfhg + -cfPeEPffc 



exp (tj/jj 



t 



Summing up for i = 1, . . . , r yields the first result. In order to prove the second statement, note that 
E[p(7?e)] < a log E [^(770)] < log a 



log I y^ exp(tjhj 



Ei B. 



< log(a). 



For each i, Uhi^g + Yi^Pi.g is a convex function of the optimization parameters 9. By [11, Example 3.14], 

given convex functions g\(9), . . . ,g r (9), the function f(9) = log(e 9l ' e ' + • ■ ■ + e 9fc ^^) is itself convex. It follows that 
log £[1/3(7/0)] is a convex function of 9 and that the constraint set 

{deQ: logE^M] < loga} = {9 EG: E[<p(r)e)} < a} 

is convex. Finally, from Lemma 9, if £[1/5(770)] < a then ¥[i]g < 0] > 1 — £[(^(7/0)] > 1 — a. Together with Assumption 1, 
this implies that the optimization problem with constraint CX P (tihi 

a convex conservative approximation of (16)— (18). 



IIIE'P, 



< a in place of (18) is 
□ 



This convex approximation of Problem (16)~(18) is obtained in (13)-(15) by setting R = 1 and (f>{rj) = <p(rj) — a. 
The result that the approximation is conservative relies essentially on the fact that, with this choice, </)(?/) > 4>cc{v) 
V?/ G M. r (see Lemma 9). This result can be generalized: Given any two functions <fi',(f>" : R r — > R such that 
ftiv) ^ 0"( r7 ) Vt/ G K r , constraint (15) with = 0' is more conservative than the same constraint with (/> = <//'. This 
type of analysis can be exploited to compare different probabilistic constraints and to minimize the conservatism of 
the convex approximations with respect with the tunable parameters, but is not fully pursued here. 

3.3 Ellipsoidal Constraint Functions 

Consider the constraint function 



n(xg,ug) 















Xg 




Xg 




( 




-)'■( 




-) 




_Ug_ 




Ug_ 





where S > and 6 are given. Then the constraint n(xg,ug) < restricts the vector xj to an ellipsoid with 

center S and shape determined by S. We now provide an approximation of the chance constraint (18) that is a 
semi-definite program in the optimization parameters 9 = (G,d). Similar to §3.2.2, the idea is to ensure that the 
100(1 — a)% confidence ellipsoid of w is such that (18) holds. To this end, let 



ye 







Xg 






-5 = 


Ug_ 





Ax + B(Gw + d) 
Gw + d 



Di 



5 = tin + PLw, 



where h' a 



Ax + Bd 
d 



S and P'n 



BG + D 
G 
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Proposition 11 Define Sg = /3(a)H2P^£s , with /9(a) as in Proposition 7, and 6 = B^h' e . Then 

-A+i o a 



XI (S e ) T 
6 Sg I 



>0, 



A > 0, 



(26) 



is a Linear Matrix Inequality (LMI) in the unknown parameters 9 and A. If (0,X) is a solution of (26), then 
satisfies (18). As a consequence, under Assumption 1, the problem 

inf E[V(x e ,ug)} 
see, a 

subject to (7), (9) and (26) 

is a convex conservative approximation of Problem (16)-(18). 

Proof: The inequality (18) may be equivalently represented as 

P(yjEy 9 - 1 < 0) > 1 - a 

P(\\~iye\\ 2 -1 <0) > I- a 



|H^^ + P»f-1<0) >1 
|6 + S>|| 2 - 1 < 0) >l-a, 



(27) 



where 5^ = H^. Since w ~ 7V(0, £), one can compute /3(a) > such that ||S 1 / 2 w|| < /3(a) 2 specifies the required 

100(1 — a) % confidence ellipsoid of w. Hence, we need to ensure that ||£~ 1/2 w|| < /3(a) 2 =>||6 + S' g w\\ 2 < 1. This 
is equivalent to 

Sup ||6 + ^W|| 2 < 1 <=> Sup ||6 + Sgvf < 1. 

||s- 1 /2 fi) || 2 < (3 ( ct )2 l|n|| 2 <i 

It follows from [11, p. 653] that supi^n^ j|6 + Sgv\\ 2 < 1 if and only if there exists A > such that 



-6 T 6-A + i $s e 

(5 e ) T 6 XI - (Se) T S s 



> 0. 



Using Schur complements the last relation can be rewritten equivalently as (26). Therefore, any solution of (26) 
implies (18). To verify that (26) is an LMI, note that Sg and 6 are affine in the optimization variables. Together 
with the assumed convexity of E[V(xe, ug)], the last statement of the proposition follows. □ 



4 Integrated chance constraints 

In this section we focus on the problem 

inf E[V(xg,ue)} (28) 

subject to (7), (9) and (29) 
Ji{e)<Pi, i = l,...,r (30) 

where, for i = 1, . . . ,r, Ji(8) = K[<pi o rji(xg,ug)}, functions cpi : K — > K are measurable and the > are fixed 
parameters. This problem corresponds to Problem (13)-(15) when setting R — r and 4>i(j]) = ipi(rf) — Pi, with 
i = 1, . . . , r. For the choice 

j0, iiz<0, 
1 \z, otherwise, 
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constraints of the form (30) are known as integrated chance constraints [19,20]. In fact, one may write (dropping the 
dependence of 771 on xe and ug to simplify the notation) 



E[<pi(r)i)]=E 



tX 



Mo.^is) ds 



(rji > s) ds, 



where ls(-) is the indicator function of set S and the second equality follows from Tonclli's theorem [16, Theorem 
4.4.5]. Therefore, constraint (30) is equivalent to 

+00 

P(r)i(x e ,u e ) > s)ds <p h (32) 

whence the name integrated chance constraint. Note that (pi plays the role of a penalty (or barrier) function that 
penalizes violations of the inequality r/i(x, u) < 0, and Pi is a maximum allowable cost in the sense of (32). Of course, 
different choices of ipi will not guarantee the equivalence between (30) and (32). However, they may be useful in 
deriving other quantitative chance constraint-type approximations. 

4-1 Convexity of Integrated Chance Constraints 

We now establish sufficient conditions on the rji and (pi for the convexity of the constraint set 

^icc = {9: Ji{0)<Ph i=l,...,r}. (33) 



The result is again a consequence of Proposition 3. 

Proposition 12 Let the mappings rji : ^t^+iJ^x^™ 1 M. be measurable and convex, and let the ipi : M — > R be 
measurable, monotone nondecreasing and convex. Assume that the Ji(9) are finite for all 9. Then each Ji(9) is 
a convex function of 9 and J^icc is a convex set. As a consequence, under Assumption 1, (28) -(30) is a convex 
optimization problem. 

Proof: Fix w G arbitrarily. Since the mapping 9 1— > {xg{uj) 1 ug{u})') is affine and (x,u) > r]i(x,u) is convex by 
assumption, their composition 9 1— > rjAxiw), u{u)j \ is a convex function of 9. Using the assumption that (pi is mono- 
tone nondecreasing and convex, we may apply Proposition 3 with 7(0;, 9) = rji{xe{oj), ug(uj)) and (fi in place of ip to 
conclude that Ji(9) = E[<pj o 7(0;, 9)] is convex. Hence, for any choice of Pi, the set ^i = {9 : Ji{9) < Pi} is convex. 
Since &icc = f~li=i the convexity of J^icc follows. Together with Assumption 1, this proves that (28)-(30) is a 
convex optimization problem. □ 



It is worth noting that the function ip(-) of Section 3.2.4 satisfies analogous monotonicity and convexity assumptions 
with respect to each of the r)i, with i = 1, . . . , r. Unlike those of Section 3, this convexity result is independent 
of the probability distribution of w. By virtue of the alternative assumptions (iii') and (hi") of Proposition 3, the 
requirement that J(9) be finite for all 9 may be relaxed. A sufficient requirement is that there exist no two values 
9 and 9' such that J{9) = +00 and J{9') = —00. In particular, provided measurable and convex rji, definition (31) 
satisfies all the requirements of Proposition 12. 

Example 13 The (scalar) polytopic constraint function: 

r)(x g ,u g ) = T x x e + T u u g - y 

fulfills the hypotheses of Proposition 12. Hence, the corresponding integrated chance constraint is convex. 

Example 14 Following Example 1, an interesting case is that of ellipsoidal constraints. For an ((N+ l)?i + A^m)-size 
positive-scmidefinite real matrix 3 and a vector 5 S ^( N + 1 ) n + Nm ; define 



rj{xe,ue) = 







T 






( 


xe 


-H 


x e 


-) 




_U0_ 




_u g _ 





(34) 
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This is a convex function of the vector 

T 



the affine mapping £ = 



Xn lla 



(it is the composition of the convex mapping £ T S£, S > 0, with 
5) and hence Proposition 12 applies. 



Remark 15 A problem setting similar to Example 14 with quadratic expected-type cost function and ellipsoidal 
constraints has been adopted in [26], where hard constraints are relaxed to expected-type constraints of the form 
£[77(5, u)\ < (3. This formulation can be seen as a special case of integrated chance constraints with <p(z) = z for all 
z. The choice of tp within a large class of functions is an extra degree of freedom provided by our framework that 
may be exploited to establish tight bounds on the probability of violating the original hard constraints, see Lemma 9 
for an example. 

4-2 Numerical solution of optimization problems with ICC 

Even though ICC problems are convex in general, deriving efficient algorithms to solve them is still a major chal- 
lenge [19,20]. For certain ICCs, it is possible however to derive explicit expressions for the gradients of the constraint 
function. Provided the cost has a simple (e.g. quadratic) form, this allows one to implement standard algorithms 
(e.g. interior point methods [11]) for the solution of the optimization problem. 

Let w satisfy Assumption 2 with, for simplicity, S equal to the identity matrix, i.e. w ~ Af(0,I). Consider the 
problem with one scalar constraint (the generalization to multiple (joint) constraints is straightforward): 



-T -T 

Xn Un 



M 



(35) 



subject to (7), (9) and E[ip(T x x B + T u u g - y)] < 0, 

where ip is defined as in (31). 

Lemma 16 Let z be a Gaussian random variable with mean \i and variance a 2 > 0. Then 



where 



^) = |erfc(-|)+- 7 =exp(^ 
and erfc(a:) = 1 — erf(cc) is the standard complementary error function. 
Proof: Since z ~ A/"(/i, cr 2 ) it holds that 



n<p(*)] 



By the change of variable y = one gets 



£exp - - h J ) dt. 



\Z2na 2 Jo 



2a 2 



E[<p(z)] = 4= / 0" + ^v) cM-y 2 )dy 

\/2cr 

= —= \ exp(-y 2 )dy + aJ- yexp(-y 2 )dy 



H c ( — At \ a 

- erfc —=- + —== exp 

2 \V2aJ V2^ 



,2 



2a 2 r 



which is equal to ag(-). 



□ 
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Simple calculations and the application of this lemma yield the following result. 

Proposition 17 Problem (35) is equivalent to 

minimize hi(d) + /12(G) 
subject to o~g( 

where 



(36) 



hi(d) 



^=(T x Ax ~y) + (T x B + T u )d, 
a=\\{T x D + {T x B + T u )G) T \\ 2 , 

Ax + Bd 



{Ax + Bd) T 3F 



M 



h 2 (G)=Tr[ 



(BG + D) T G T 



M 



d 

BG + D 
G 



Note that expectations have been integrated out. Now it is possible to_put the problem in a standard form for 
numerical optimization. Let Gi and Di be the i-th column of G and D, respectively. Redefine the optimization 



variable 6 as the vector 6 



d T G\ ■ ■ ■ G% r 



. Define e, /, v and L by 



fi = (T x Ax -y) + (T X B + T u )d =e + f T 6 
(T X D + (T X B + T U )G) T =v + L0. 



Corollary 18 Problem (35) is equivalent to 



minimize 



fo(0) =3o(0)+$>(0) 



subject to H9 =0, 
h{0) <0 

where the matrix H in the equality constraint accounts for the causal structure of G, while 



90(0) =( AX J 
9i(0) 

f 1 (6)=\\v + L6 





~B~ 


T 


Ax 




~B~ 


+ 




H 9j M( 


+ 






I 







I 



H e), 







~B~ 





+ 


I 



T 

Hi 9) M 



e + f 



D, 


/?<0, 



Ha 



\v + L6\\ 

and Hq and Hi are selection matrices such that Hq9 = d and HiO = Gi 



We conclude the section by documenting the expressions of the gradient and the Hessian of the constraint function 

V7*(5\ 1 tT v + L8 ( (e + f T 9f \ 1 / e + f T S \ 
V/i(«) = — £ ii — — ^expf —^) +-erfc( ^ — — r ^ i )f 



^2tt \\v + L0\ 



2\\v + L6\\ 



■ crfc ( 

2 V y/2 \\v + L9\ 
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Fig. 2. A mechanical system with springs and masses. 



* AW _ i (_i£±£51 

V2^ V 2\\v + L6\ 



+ MS)- MS) 



where 



Ji(0) = 
MS) = 
MS) = 



i 



\v + L6 



-{L T L + ff T ), 



(e + / J 6») 2 - \\v + L6\ 



L0\ 



e + f T Q 
\v + L6\f 



(L T (v + L8){v + L9) T L), 



(L T (v + L9)f T + f(v + L6) T L). 



The expressions of gradient and Hessian of the quadratic function fo(6), used e.g. by interior point method solvers, 
are quite standard and will not be reported here. 



5 Simulation results 



We illustrate some of our results with the help of a simple example. Consider the mechanical system shown in 
Figure 2. di, • • • , d± are displacements from an equilibrium position, u\, ■ ■ ■ ,1*3 are forces acting on the masses. In 
particular, u\ is a tension between the first and the second mass, u-i is a tension between the third and the fourth 
mass, and 1*3 is a force between the wall (at left) and the second mass. We assume all mass and stiffness constants 
to be equal to unity, i.e. mi = • • ■ = TO4 = 1, fci = • ■ ■ = fei = 1. We consider a discrete-time model of this system 
with noise in the dynamics, 



x(t + 1) = Ax(t) + Bu(t) + w(t), 



where w is an i.i.d. noise process, w(t) ~ J\f(0,<j 2 I),<r = 0.05 for all t, and x = d\, d®, (I3, di, di, tfe, cfe, ^4 ■ The 

discrete-time dynamics arc obtained by uniform sampling of a continuous-time model at times t ■ h, with sampling 
time h = 1 and t = 0, 1,2, . . ., under the assumption that the control action u(t) is piecewise constant over the 
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sampling intervals [t ■ h, t ■ h + h). Hence, A = e c and B = A c 1 (e hAc - — T)B C , where A c and B c , denned as 



A n 



04x4 Ax4 
A C 21 04x4 



04x3 
B c ,2l 



A Ct2 l = 



B r 



-2 1 0, 
1-210 
1-21 

1-1 

1 
-1 -1 
1 
0-10 



are the state and input matrices of the standard ODE model of the system. 

We are interested in computing the control policy that minimizes the cost function 



E 



N-l 



(x{t) T Qx(t) + u(t) T Ru(t)) + x{N) T Qx(N) 



t=o 



where the horizon length is fixed to N — 5, the weight matrices are defined as Q 



I 




(penalizing displacements 



but not their derivatives) and R = I. The initial state is set to xq = [0, 0, 0, 1, 0, 0, 0, 0] r . In the absence of constraints, 
this is a finite horizon LQG problem whose optimal solution is the linear time-varying feedback from the state 

u(t) = -(B T P(t + l)B + R)~ 1 B T P(t + l)Ax(t), 

where the matrices P(t) are computed by solving, for t = N— 1, . . . , 0, the backward dynamic programming recursion 

P(t) = Q + A T P(t + I) A - A T P{t + l)B{B T P{t + l)B + R)- 1 B T P{t + 1)A, 

with P(N) = Q. Simulated runs of the controlled system are shown in Figure 3. We shall now introduce constraints 
on the state and the control input and study the feasibility of the problem with the methods of Section 3. The 
convex approximations to the chance-constrained optimization problems are solved numerically in Matlab by the 
toolbox CVX [18]. In all cases we shall compute a 5-stage affine optimal control policy and apply it to repeated runs 
of the system. Based on this we will discuss the feasibility of the hard constrained problem and the probability of 
constraint violation. 



5.1 Polytopic constraints 



Let us impose bounds on the control inputs, < 0.1, |«2(i)| < 0.3 and |i*3(t)| < 0.15, with t = 0, . . . , N — 1, 

and bounds on the mass displacements, |di(f)| < 10, for i = 1, • • • , 4 and with t = 1, . . . , N. In the notation of 

1 T 



Section 3.2, these constraints are captured by the equation rj(x, u) = T x x + T u u — y < where T x 

1 T 



M 1 



and 



H 1 



with 



M 



Ah 



Mi 



^4x4 04x4 




' I ' 






_~ ^4x4 04x4_ 




-I 
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time 



Fig. 3. 1000 sample paths of the system with LQG control policy. Above: control input. Below: mass displacements. 



and 



Vi = 



10 




V 




"o.io 




, V2 = 




, y' = 


0.30 


10 




y'_ 




0.15 



This hard constraint is relaxed to the probabilistic constraint ¥[rj(x,u) < 0] > 1 — a. The resulting optimal control 
problem is then addressed by constraint separation (Section 3.2.1) and ellipsoidal approximation (Section 3.2.2). 

With constraint separation, the problem is feasible for a > 0.05. For a = 0.1, the application of the suboptimal 
control policy computed as in Proposition 6 yields the results shown in Figure 4. With this policy, the control input 
saturates within the required bounds whereas the mass displacements stay well below bounds. In fact, although the 
required probability of constraint satisfaction is 0.9, constraints were never violated in 1000 simulation runs. This 
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Fig. 4. 1000 sample paths of the system with control policy computed via constraint separation. Above: control input. Below: 
mass displacements. Horizontal straight lines show bounds. 



suggests that the approximation incurred by constraint separation is quite conservative, mainly due to the relatively 
large number of constraints. It may also be noticed that the variability of the applied control input is rather small. 
This hints that the computed control policy is essentially open-loop, i.e. the linear feedback gain is small compared 
to the affine control term. 

With the ellipsoidal approximation method, for the same probability level, the problem turns out to be infcasiblc, 
in accordance with the conclusions of Section 3.2.3. For the sake of investigation, we loosened the bounds on the 
mass displacements to < 100 for all i and t. The problem of Proposition 7 is then feasible and the results 

from simulation of the controlled system are reported in Figure 5. Although the controller has been computed under 
much looser bounds, the control performance is similar to the one obtained with constraint separation, a clear sign 
that the ellipsoidal approximation is overly conservative in this case. Another evidence of inaccuracy is the fact that, 
while the control inputs get closer to the bounds, the magnitude of the displacements is not reduced. As in the case 
of constraint separation, the applied control input is insensitive to the specific simulation run, i.e. the control policy 
is essentially open loop. 

5.2 Ellipsoidal constraints 

Consider the constraint function 



with c = (l 2 + l 2 + l 2 + l 2 + 0.1 2 + 0.15 2 + 0.3 2 ) = 4.1225. Unlike the previous section, we do not impose bounds on 
di{k) and m(k) at each k but instead require that the total "spending" on x and u does not exceed a total "budget". 
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Fig. 5. 1000 sample paths of the system with control policy computed via ellipsoidal approximation. Above: control input. 
Below: mass displacements. Horizontal straight lines show bounds. 

This constraint can be modelled in the form of Section 3.3, namely 

ri(x,u) = 




1 <0, 



with 5 



and S = -j^— 

N-c 



L 
I 



where 



L 



Li 



Li 



^4x4 04x4 
04x4 04x4 



The constrained control policy for N = 5 and a = 0.1 is computed by solving the LMI problem of Proposition 11. 
Results from simulations of the closed-loop system are reported in Figure 6. Once again, constraints were not 
violated over 1000 simulated runs, showing the conservatism of the approximation. It is interesting to note that 
the displacements of the masses are generally smaller than those obtained by the controller computed under affinc 
constraints, at the cost of a slightly more expensive control action. In contrast with the affine constraints case, the 
control action obtained here is much more sensitive to the noise in the dynamics, i.e. the feedback action is more 
pronounced. 
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Olr- 




0.5 1 1.5 3 2.5 3 3.5 4 4.5 5 



Fig. 6. 1000 sample paths of the system with control policy computed via ellipsoidal constraints. Above: control input. Below: 
mass displacements. Horizontal straight lines show bounds. 

6 Conclusions 

We have studied the convexity of optimization problems with probabilistic constraints arising in model predictive 
control of stochastic dynamical systems. We have given conditions for the convexity of expectation-type objective 
functions and constraints. Convex approximations have been derived for nonconvex probabilistic constraints. Results 
have been exemplified by a numerical simulation study. 

Open issues that will be addressed in the future are the role of the tunable parameters (e.g. the cti in Section 3.2.1, 
the ti of Section 3.2.4 and the /3, in Section 4) in the various optimization problems, and the effect of different choices 
of the ICC functions ifi (Section 4). Directions of future research also include the extension of the results presented 
here to the case of noisy state measurements, the exact or approximate solution of the stochastic optimization 
problems in terms of explicit control laws and the control of stochastic systems with probabilistic constraints on the 
state via bounded control laws. 
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